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A  general  electrode-electrolyte-assembly  (EEA)  model  has  been  developed,  which  is  valid  for  different 
designs  of  solid  oxide  fuel  cells  (SOFCs)  operating  at  different  temperatures.  In  this  study,  it  is  applied  to 
analyze  the  performance  characteristics  of  planar  anode-supported  SOFCs.  One  of  the  novel  features  of 
the  present  model  is  its  treatment  of  electrodes.  An  electrode  in  the  present  model  is  composed  of  two 
distinct  layers  referred  to  as  the  backing  layer  and  the  reaction  zone  layer.  The  other  important  feature  of 
the  present  model  is  its  flexibility  in  fuel,  having  taking  into  account  the  reforming  and  water-gas  shift 
reactions  in  the  anode.  The  coupled  governing  equations  of  species,  charge  and  energy  along  with  the 
constitutive  equations  in  different  layers  of  the  cell  are  solved  using  finite  volume  method.  The  model  can 
predict  all  forms  of  overpotentials  and  the  predicted  concentration  overpotential  is  validated  with  mea¬ 
sured  data  available  in  literature.  It  is  found  that  in  an  anode-supported  SOFC,  the  cathode  overpotential 
is  still  the  largest  cell  potential  loss  mechanism,  followed  by  the  anode  overpotential  at  low  current  den¬ 
sities;  however,  the  anode  overpotential  becomes  dominant  at  high  current  densities.  The  cathode  and 
electrolyte  overpotentials  are  not  negligible  even  though  their  thicknesses  are  negligible  relative  to  the 
anode  thickness.  Even  at  low  fuel  utilizations,  the  anode  concentration  overpotential  becomes  significant 
when  chemical  reactions  (reforming  and  water-gas  shift)  in  the  anode  are  not  considered.  A  parametric 
study  has  also  been  carried  out  to  examine  the  effect  of  various  key  operating  and  design  parameters  on 
the  performance  of  an  anode-supported  planar  SOFCs. 

©  2008  Elsevier  B.V.  All  rights  reserved. 


1.  Introduction 

Solid  oxide  fuel  cells  (SOFCs)  are  receiving  considerable 
attention  for  stationary  power  applications  due  to  their  high 
energy  efficiency,  low  pollutant  emissions,  fuel  flexibility,  and 
co-generation  capability.  Among  various  configurations  in  which 
SOFCs  exist,  two  common  configurations  are  tubular  and  planar. 
Although  significant  progresses  have  been  made  for  the  tubular 
configuration  of  SOFC,  tubular  SOFC  has  low  power  density  due 
to  high  electrical  resistance  arising  from  long  current  paths  [1].  In 
contrast,  the  planar  configuration  of  SOFC  is  capable  of  achieving 
very  high  power  density  [1,2].  However,  due  to  high  temperature 
operation,  there  exist  technical  challenges  for  planar  SOFCs,  such  as 
requirement  of  high  temperature  gas  seals,  high  material  and  man¬ 
ufacturing  costs.  Therefore,  it  is  important  to  reduce  the  operating 
temperature  to  an  intermediate  range  (823-1073  K),  often  referred 
to  as  intermediate-temperature  SOFCs;  then,  conventional  stainless 
steel  can  be  used  for  interconnects  instead  of  more  expensive  high 
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chrome  alloys  or  oxides,  and  the  material  and  manufacturing  costs 
can  be  reduced  significantly  [3,4]. 

Planar  SOFCs  are  generally  fabricated  into  two  basic  designs, 
namely  electrolyte-supported  and  electrode-supported.  In 
electrolyte-supported  SOFCs,  electrolyte  is  the  thickest  com¬ 
ponent  on  which  anode  and  cathode  layers  are  deposited.  Since 
electrolyte  is  the  thickest  component,  the  ohmic  contribution  of 
the  electrolyte  layer  is  significant  in  electrolyte-supported  SOFCs; 
therefore,  they  are  applicable  for  high  temperature  operation 
with  the  typical  electrolyte  material,  which  is  yttria  stabilized 
zirconia  (YSZ)  [5].  On  the  other  hand,  in  electrode-supported 
SOFCs,  either  anode  or  cathode  is  the  thickest  component  on  which 
all  other  layers  are  deposited.  Among  electrode-supported  SOFCs, 
anode-supported  design  is  preferred  over  cathode-supported 
design  for  operation  at  reduced  temperatures.  This  is  due  to  higher 
contribution  of  cathode  overpotentials  to  the  total  cell  potential 
loss  than  the  anode,  especially  at  reduced  operating  temperatures 
[5,6]. 

The  processes  influencing  an  SOFC  performance  in  different  lay¬ 
ers  of  the  cell  are  complex,  competing  and  interdependent.  Because 
of  expensive  and  time  consuming  physical  prototyping,  experimen¬ 
tal  advances  are  quite  limited.  Hence,  to  enhance  the  development 
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Nomenclature 

Av  reactive  surface  area  per  unit  volume  (m2  m-3 ) 

Ch2  hydrogen  concentration  (mole  m-3 ) 

Ch2  ref  reference  hydrogen  concentration  (mole  m3 ) 

Co2  oxygen  concentration  (mole  m-3 ) 

c0 2  ref  reference  oxygen  concentration  (mole  m-3 ) 
cpi  specific  heat  of  species  i  at  constant  pressure 
(J  mole-1 1<-1 ) 

dp  diameter  of  the  pore  ( pun) 

Dy  binary  diffusion  coefficient  (m2  s-1 ) 

DKn  i  Knudsen  diffusion  coefficient  (m2  s-1 ) 
effective  diffusion  coefficient  (m2  s-1 ) 

F  Faraday’s  constant  (96,487  C  mole-1 ) 

hf  absolute  enthalpy  of  species  i  (J  mole-1 ) 

J  current  density  (A  m-2 ) 

Je  electronic  current  density  (A  m-2 ) 

Ji  ionic  current  density  (A  m-2 ) 

Jo2ref  reference  exchange  current  density  for  H2  oxidation 

(Am-2) 

Jo2ref  reference  exchange  current  density  for  02  reduction 

(Am-2) 

k  thermal  conductivity  (W  m-1 1<-1 ) 

kei  thermal  conductivity  of  pure  electron-conducting 
particles  (W  m-1  K-1 ) 

kf  thermal  conductivity  of  fluid  mixture  (W  m-1  K-1 ) 

kf  thermal  conductivity  of  pure  component  i 

(Wm-1 1<-1) 

ki0  thermal  conductivity  of  pure  electron-conducting 
particles  (W  m-1  K-1 ) 

ks  thermal  conductivity  of  solid  phase  of  backing  layers 

(Wm-1 1<-1) 

kbr  backward  reaction  rate  constant  for  reforming  reac¬ 
tion  (mole  m-3  Pa-2  s-1 ) 

kfr  forward  reaction  rate  constant  for  reforming  reac¬ 
tion  (mole  m-3  Pa-2  s-1 ) 

kbs  backward  reaction  rate  constant  for  shift  reaction 
(molem-3  Pa-2  s-1) 

kfS  forward  reaction  rate  constant  for  shift  reaction 
(molem-3  Pa-2  s-1) 

I<p  equilibrium  constant  for  partial  pressure 
Kpr  equilibrium  constant  for  reforming  reaction 
Kps  equilibrium  constant  for  shift  reaction 
Mi  molecular  weight  of  species  i  (kg  mole-1 ) 
n  moles  of  electrons  transferred  per  mole  reactant 

nei  number  fraction  of  electron-conducting  particles  in 

the  reaction  zone  layers 

ni0  number  fraction  of  ion-conducting  particles  in  the 
reaction  zone  layers 

nt  total  number  of  particles  in  the  reaction  zone  layers 

per  unit  volume  m-3) 

Nj  molar  flux  of  species  i  (mole  m-2  s-1 ) 

p  pressure  (Pa) 

pei  probability  of  percolation  of  electron-conducting 

particles  in  reaction  zone  layers 
pio  probability  of  percolation  of  ion-conducting  parti¬ 

cles  in  reaction  zone  layers 
Pi  partial  pressure  of  species  i  (Pa) 

q  heat  flux  ( W  m-2 ) 

rel  radius  of  electron-conducting  particles  in  the  reac¬ 

tion  zone  layers  (m) 

ri0  radius  of  ion-conducting  particles  in  the  reaction 

zone  layers  (m) 


r{  rate  of  reaction  (mole  m-2  s-1 ) 

rr  volumetric  reforming  reaction  rate  (kg  m-3  s-1 ) 

rs  volumetric  shift  reaction  rate  (kg  m-3  s-1 ) 

R  universal  gas  constant  (8.3143  J  mole-1 1<-1 ) 

JZa  volumetric  current  density  produced  due  to  hydro¬ 
gen  oxidation  (Am-3) 

1ZC  volumetric  current  density  produced  due  to  oxygen 
reduction  (Am-3) 

Ash2  molar  entropy  change  for  an  overall  fuel  cell  reaction 

involving  hydrogen  (J  mole-1 1<-1 ) 

As0  molar  entropy  change  for  a  reaction  at  standard 
pressure  (J  mole-1  K-1 ) 

Se  rate  of  production  or  consumption  of  energy 

(Wm-3) 

Ss  i  rate  of  production  or  consumption  of  species  i 
(molem-3  s-1) 

T  temperature  (I<) 

x,  mole  fraction  of  species  i 

Z  average  coordination  number 

Zel  coordination  number  of  electron-conducting  parti¬ 

cles  in  the  reaction  zone  layers 
ZI0  coordination  number  of  ion-conducting  particles  in 

the  reaction  zone  layers 

Greek  symbols 

a  charge  transfer  coefficient 

Yh2  reaction  order  for  hydrogen  oxidation  reaction 

e  porosity 

p  electrode  activation  overpotential  (anode  or  cath¬ 

ode)  (V) 

0  contact  angle  between  electron-  and  ion¬ 

conducting  particles  in  the  reaction  zone  layers 
k  ionic  conductivity  (S  m-1 ) 

/x,  dynamic  viscosity  of  pure  component  i  (kg  m-1  s-1 ) 

vh2  stoichiometric  coefficient  of  H2  involved  in  the  oxi¬ 

dation  reaction 

vo2  stoichiometric  coefficient  of  OH2  involved  in  the 

reduction  reaction 

£  fraction  of  entropy  change  of  the  overall  SOFC  reac¬ 

tion 

a  electronic  conductivity  (S  m-1 ) 

r  tortuosity 

4>e  electronic  potential  (V) 

0,  ionic  potential  (V) 

0  volume  fraction  of  electron-conducting  particles  in 

the  reaction  zone  layers 
@ij  dimensionless  function 

Subscripts 

bl  backing  layer 

el  electrolyte  layer 

ref  reference 

rl  reaction  zone  layer 

Superscripts 
eff  effective 


and  understanding  of  SOFCs,  it  is  imperative  to  develop  a  math¬ 
ematical  model  which  accounts  for  complex  transport  processes 
with  chemical  and  electrochemical  reactions.  The  present  model 
can  be  used  to  predict  cell  performance  and  gain  insights  on  the 
physical,  chemical  and  electrochemical  processes  occurring  in  dif¬ 
ferent  layers  of  the  cell. 
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Over  the  past  few  years,  many  researchers  have  developed 
anode-supported  SOFC  models  [3,5,7-11],  which  can  be  broadly 
classified  into  macro-  or  micro-models,  depending  on  the  treat¬ 
ment  of  electrodes.  In  macro-models,  electrodes  are  modeled 
as  porous  structures  of  purely  electron-conducting  particles  and 
electrochemical  reactions  are  considered  to  occur  exclusively  at 
the  electrode/electrolyte  interfaces,  thereby  treating  reaction  zone 
layers  as  interfacial  boundary  conditions  [3,5,8-11].  Whereas,  in 
micro-models,  electrodes  are  modeled  as  porous  structures  of 
electron-  and  ion-conducting  particles.  Essentially,  electrodes  in 
micro-models  are  treated  as  reaction  zone  layers  having  triple 
phase  boundaries  (TPBs)  scattered  throughout  the  electrodes 
[7,12-15,41].  However,  most  of  the  micro-models  are  electrode¬ 
level  models  predicting  the  electrochemical  characteristics  of  an 
electrode,  either  anode  or  cathode;  only  recently  Nam  and  Jeon 
[7]  integrated  this  approach  into  a  cell-level  model.  Additionally, 
Nam  and  Jeon  [7]  and  many  others  [12,15-18]  concluded  that  the 
TPBs  are  most  active  at  the  electrode/electrolyte  interface  and  most 
of  the  electrochemical  reactions  in  the  electrodes  occur  within  a 
distance  of  the  order  of  10  p,m  from  the  electrolyte.  Conversely, 
most  of  the  macro-models  are  cell-level  models  predicting  the  cell 
performance  at  different  operating  and  design  conditions.  Incorpo¬ 
ration  of  micro-characteristics  of  electrodes  into  macro-models  not 
only  helps  in  better  understanding  the  relevant  processes  occurring 
within  the  electrodes,  but  also  enhances  the  predicting  capabil¬ 
ity  of  the  overall  cell  model.  One  of  the  novelties  of  the  present 
model  is  its  treatment  of  the  electrodes.  An  electrode  in  the  present 
model  is  treated  as  two  distinct  layers  referred  to  as  the  backing 
layer  and  the  reaction  zone  layer  for  electrochemical  reactions,  thus 
serving  as  a  bridge  connecting  the  micro-  and  macro-approach  of 
modeling  electrodes.  The  other  important  feature  of  the  present 
model  is  its  flexibility  in  fuel  choice;  not  only  pure  H2  but  also  any 
reformate  composition  (H2,  H20,  CO,  C02  and  CH4)  can  be  used 
as  a  fuel.  Further,  the  modified  Stefan-Maxwell  equations  incor¬ 
porating  Knudsen  diffusion  are  used  to  model  multi-component 
diffusion  in  the  porous  electrodes.  Moreover,  chemical  reactions 
such  as  water-gas  shift  reaction  and  methane  reforming  reaction 
are  considered  in  the  anode  of  an  SOFC. 

The  objective  of  this  study  is  to  develop  a  mathematical  model 
of  an  SOFC  that  combines  the  micro-  and  macro-models  and  takes 
into  account  all  the  relevant  physical,  chemical  and  electrochemi¬ 
cal  processes  in  different  layers  of  the  cell  and  predicts  all  forms  of 
overpotentials.  The  developed  model  is  then  applied  to  analyze  per¬ 
formance  characteristics  of  planar  anode-supported  SOFCs.  Finally, 
the  effect  of  various  key  operating  and  design  conditions  on  the 
performance  of  an  anode-supported  SOFC  is  investigated. 

2.  Model  formulation 

A  schematic  of  a  planar  anode-supported  SOFC  is  shown  in  Fig.  1 . 
Since  the  geometry  is  symmetric  about  the  middle-axis  of  the  flow 
channel,  the  schematic  shown  is  the  right  half  cross-section  of  the 
cell.  A  typical  anode-supported  SOFC  consists  of  a  thick  anode,  a 
relatively  thinner  electrolyte  and  cathode.  However,  in  the  present 
formulation,  two  finite  distinct  layers  are  considered  between  the 
electrodes  and  the  electrolyte,  which  are  referred  to  as  the  anode 
reaction  zone  layer  and  the  cathode  reaction  zone  layer,  respec¬ 
tively.  The  computational  domain  is  shown  by  a  dashed  line,  which 
includes  the  land  portion  of  the  interconnects  interfacing  the  back¬ 
ing  layers,  the  porous  portion  of  the  backing  layers  interfacing  the 
flow  channels,  the  backing  layers,  the  reaction  zone  layers,  and  the 
electrolyte  layer. 

The  cell  is  assumed  to  operate  under  steady  state  and  the  param¬ 
eters  vary  in  the  x  (across  the  cell)  and  y  (transverse)  directions. 
The  convective  flux  is  negligible  in  the  porous  backing  and  reac¬ 
tion  zone  layers  when  compared  to  the  diffusive  flux  of  gaseous 
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Fig.  1.  Schematic  of  an  anode-supported  planar  solid  oxide  fuel  cell  considered  in 
this  study. 


species,  which  means  the  primary  mode  of  species  transport  in 
the  porous  backing  and  reaction  zone  layers  is  by  diffusion  [5].  For 
instance,  according  to  Yakabe  [4],  the  calculated  diffusion  flux  of  H2 
at  1  A/cm2is  about  0.23  mole/m2s,  which  is  four  orders  of  magni¬ 
tude  higher  than  convective  flux  (1  x  10-4  mole/m2s).  The  reactant 
gas  mixtures  are  approximated  as  ideal  gases  with  negligible  Soret, 
Dofour  and  gravity  effects.  Since  the  reaction  zone  layers  are  con¬ 
sidered  as  separate  regions,  there  are  no  electrochemical  reactions 
(either  oxidation  or  reduction)  in  the  backing  layers.  Moreover,  in 
the  present  study,  it  is  assumed  that  there  is  no  electrochemical 
oxidation  of  CO  in  the  anode  reaction  zone  layer.  In  addition,  the 
electrolyte  layer  is  assumed  to  be  a  dense  solid  such  that  reactant 
gas  crossover  can  be  neglected.  With  these  assumptions,  the  cell 
model  is  formulated  and  described  in  the  following  sections. 

2.1.  Backing  Layers 

Backing  layers  are  porous  structures,  which  provide  flow  paths 
for  the  reactants  to  and  products  from  the  reaction  zone  layers. 
They  also  conduct  electrons  through  their  solid  portion,  thereby 
providing  a  flow  path  for  the  transport  of  electrons.  Therefore, 
the  processes  that  need  to  be  modeled  in  the  backing  layers  are 
transport  of  multi-component  mixture  of  species  to  and  from  the 
reaction  zone  layers  along  with  chemical  reactions  (water-gas  shift 
and  methane  reforming)  in  the  anode  backing  layer,  transport  of 
electrons  in  the  solid  portion  of  the  porous  structures,  and  trans¬ 
fer  of  energy  due  to  heat  conduction,  species  diffusion  and  heat 
generation/consumption  due  to  chemical  reactions  and  electron 
migration. 

The  mathematical  model  governing  the  processes  in  the  backing 
layers  is  formulated  by  applying  conservation  of  species  equation 
along  with  the  modified  Stefan-Maxwell  equations  incorporating 
Knudsen  diffusion  for  multi-component  gas  diffusion,  conserva¬ 
tion  of  electronic  charge,  and  conservation  of  energy.  The  governing 
equations  are  described  as  follows: 


Species  :  V  •  7?,-  =  Ss  { 

(i) 

Electronic  charge  :  V  •  ~f  e  =  0 

(2) 

Energy  :  V  -  if  =  Se 

(3) 
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where  i  =  CH4,  H2,  H20,C0,  and  C02  for  anode;  and  i  =  02,  and  N2 
for  cathode,  7?  j  is  the  molar  diffusion  flux  of  species  i,  Ss  i  is  the 
species  source  term,  ~fe  is  the  electronic  current  density,  if  is  the 
energy  flux,  and  Se  is  the  energy  source  term. 

The  diffusion  flux  in  the  porous  layers  can  be  determined  using 
the  modified  Stefan-Maxwell  equations  incorporating  Knudsen  dif¬ 
fusion  for  multi-component  systems  involving  n  species  [19]: 


n 


where  c  is  the  concentration  of  the  mixture,  xz  is  the  mole  fraction  of 
species  i,  ff  f  is  the  diffusion  flux  of  species  i,  and  £>jjff  is  the  effective 
diffusion  coefficient,  defined  as 


£>eff  _  £  (  DijDI<n,i  \ 

V  X  y  Dtj  +  DKn  i  J 


(5) 


where  Dy  is  the  ordinary  or  binary  diffusion  coefficient  and  DKni  is 
the  Knudsen  diffusion  coefficient  of  species  i,  s  and  r  are  porosity 
and  tortuosity  of  the  porous  layers,  respectively. 

Since  there  is  no  chemical  reaction  in  the  cathode  backing  layer, 
the  species  source  term  in  Eq.  (1)  corresponds  to  zero;  whereas,  in 
the  anode  backing  layer,  the  species  source  term  represents  the  rate 
of  production  or  consumption  of  species  due  to  methane  reforming 
and  water-gas  shift  reactions,  as  given  below: 

Methane  reforming  reaction: 

CH4  +  H20  ^  3H2  +  CO  AH298  =  206  kj  mole-1  (6) 


Water-gas  shift  reaction: 

C0  +  H20^H2+C02  AH298  = -41.1  kj  mole-1  (7) 


The  volumetric  reaction  rates  for  the  above  reactions  can  be 
written  as 


rr  =  k/rPCH4PH20  -  kbrPH2PC0  (8) 

i's  =  k/sPcoPH2o  -  hsPC02PH2  (9) 

where  rr  and  rs  are  volumetric  reaction  rates  for  the  methane 
reforming  and  water-gas  shift  reactions  expressed  in  mole  m-3  s-1 , 
kfr  and  kbr  are  the  forward  and  backward  reaction  rate  constants  for 
the  reforming  reaction  and  kfS  and  kbs  are  the  forward  and  backward 
reaction  rate  constants  for  the  water-gas  shift  reaction,  respec¬ 
tively.  In  terms  of  mole  fractions,  Eqs.  (8)  and  (9)  can  be  written 
as 


rr  =  P2lkfrXCH4XH2o  -P2kbrxl2xCo] 

(10) 

rs  =  P2[kfsxc 0Xh2o  -  kbsxc q2xH2  ] 

(ID 

where  p  is  the  total  pressure.  The  equilibrium  constants  for  the 
methane  reforming  and  water-gas  shift  reactions  can  be  deter¬ 
mined  from  the  following  empirical  relations  [20]: 

Kpr  =  =  1.0267  x  1010  x  exp(-0.2513£4  +  0.3665f3 

kbr 

+  0.5810f2  -27.134^  +  3.2770)  (12) 


Kps  =  ^=  exp(— 0.2935f3  +  0.635C2  +  4.1788f  +  0.3169) 
kbs 

where 


(= 


1000 

WJ 


(13) 


(14) 


The  forward  reaction  rate  constants  for  the  methane  reforming  and 
water-gas  shift  reactions  are  given  as  [20]: 


kfi-  =  2395  exp  (-  23^66  )  (mole  m-3  Pa-2  s-1)  (15) 


kfS  =  0.0171  exp 


103,191  J  mole-1  \ 

RT  J 


(molem  3  Pa  2  s  1 ) 


(16) 


where  R  is  the  universal  gas  constant  and  T  is  the  temperature. 

Knowing  the  forward  reaction  rate  constants  and  equilibrium 
constants  for  the  methane  reforming  and  water-gas  shift  reactions, 
the  backward  reaction  rate  constants  can  be  obtained  using  the  first 
equality  in  Eqs.  (12)  and  (13). 

Finally,  the  rates  of  production  or  consumption  of  various  species 
in  the  anode  backing  layer  can  be  formulated  as  follows: 


5s,ch4  =  -rr  (17) 

i,H2  —  3rr  +  rs  (18) 

5s,h2o  =  -rr  -  rs  (19) 

Ss,c  o  =rr-rs  (20) 

5s,co2  =  rs  (21) 

Using  Ohm’s  law,  the  electronic  current  density  given  in  Eq.  (2) 
can  be  written  in  terms  of  electronic  potential,  then  the  conserva¬ 
tion  of  charge  in  the  backing  layers  becomes 

Electronic  charge  :  V  •  (cry[fV0e)  =  0  (22) 

where  <r^ff  is  the  effective  electronic  conductivity  of  the  backing 
layers,  defined  as 

ab?  =  (f~r) a  C23) 

where  £  and  r  are  the  porosity  and  tortuosity  of  the  porous  backing 


layers,  respectively,  and  a  is  the  electronic  conductivity  of  the  pure 
backing  layer  material  (or  the  bulk  electronic  conductivity). 

The  heat  flux  “if  in  Eq.  (3)  accounts  for  the  contributions  result¬ 
ing  from  interdiffusion  of  multi-component  mixture  in  addition  to 
the  conductive  flux,  which  is  expressed  as 

n 

If  = -kl?VT +  Y)hrt  (24) 

1  =  1 

where  hz  is  the  absolute  enthalpy  of  the  species  i,  and  k%f  is  the 
effective  thermal  conductivity  of  the  backing  layers,  defined  as 

kbif  =  ek/  +  n  ~e)ks  (25) 

where  kf  is  the  thermal  conductivity  of  fluid  mixture  in  the  back¬ 
ing  layers  and  ks  is  the  thermal  conductivity  of  the  solid  phase  of 
the  porous  backing  layers.  The  thermal  conductivity  of  the  multi- 
component  fluid  mixture  is  obtained  as  [21  ] 


(26) 


where  xz  is  the  mole  fraction  of  species  i  in  a  mixture,  /q  is  the 
thermal  conductivity  of  the  pure  component  i,  and  0^  is  defined  as 


0ii  V8 


-(1/2) 


1  + 


(27) 


where  Mz  is  the  molecular  weight  of  species  i,  and  /xz  is  the  dynamic 
viscosity  of  pure  component  i. 

The  energy  source  term  in  Eq.  (3)  accounts  for  Joule  heating 
effect  representing  the  resistance  to  electron  flow  in  the  solid  phase 
of  porous  backing  layers  and  energy  generation  or  consumption 
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Table  1 

Governing  equations  in  the  backing  layers. 


Governing  equations 

Anode  backing  layer 

Cathode  backing  layer 

Species 

V^,  =  ss>i 

V  •  T$i  =  0 

Modified  S-M 

vxi = e;, 

Electronic  charge 

V  •  =  0 

v  •  (<V0O  =  0 

Energy 

v  •  (tfvr)  -  {ElMJ  ■  VT+ |f  -  (ELM  O  = 0 

v  •  (ka vr)  -  (ELicp^)  • VT + ;§ = 0 

due  to  chemical  reactions  in  the  anode  backing  layer,  which  can  be 
expressed  as 


Se  = 


j,2 


e 

eff 


(28) 


Substituting  the  expressions  of  energy  flux  and  energy  source  term 
in  Eq.  (3),  the  final  form  of  energy  equation  in  the  backing  layers 
becomes 


v-(k^vr) 


f  2 

.  VT  +  — - 

creff 

°bl 


=  0 


(29) 


where  cpi  is  the  specific  heat  at  constant  pressure  for  the  species  i. 

It  is  worthwhile  to  mention  that  the  last  term  on  the  left-hand 
side  of  above  equation  vanishes  for  cathode  backing  layer  since 
there  is  no  chemical  reaction  in  the  cathode.  The  governing  equa¬ 
tions  in  the  backing  layers  are  summarized  in  Table  1 . 


2.2.  Reaction  zone  layers 

Reaction  zone  layers  are  relatively  thin  layers  having  thicknesses 
of  the  order  of  10  [xm.  These  are  the  regions  where  fuel  and  oxi¬ 
dant  are  electrochemically  converted  into  electrical  work,  heat 
and  water  vapor.  Since  both  electron-  and  ion-conducting  parti¬ 
cles  co-exist  with  reactant  species,  these  regions  can  be  treated  as 
composite  electrodes. 

Reaction  zone  layers  have  been  modeled  as  thin  film  models, 
random  resistor  network  models,  random  packing  sphere  mod¬ 
els  or  macroscopic  porous-electrode  models,  respectively.  The  thin 
film  models  [22,23]  are  based  on  the  assumption  that  the  three 
phases  of  reaction  zone  layers  (i.e.,  electronic,  ionic,  and  gas) 
form  straight  paths  from  the  electrolyte  to  the  backing  layers. 
Because  of  this  assumption  in  thin  film  models,  the  reaction  zone 
layers  are  represented  by  a  very  ordered  structure  in  place  of 
highly  disordered  structure  as  observed  in  the  experimental  fea¬ 
tures  of  composite  electrodes  [15].  In  random  resistor  network 
models  [24,25],  the  composite  electrode  is  assumed  to  consists 
of  electron-  and  ion-conducting  particles  packed  together  to  form 
a  continuous  network.  These  models  involve  storage  of  informa¬ 
tion  in  terms  of  individual  particle  locations  and  connectivity 
with  other  particles  to  construct  resistor  networks,  followed  with 
calculation  of  electric  potential  at  each  particle  location.  Due 
to  high  computational  costs,  the  use  of  random  resistor  net¬ 
work  models  is  limited  to  problems  with  small  specimen  size 
[7].  In  random  packing  sphere  models  [13,14,41],  the  composite 
electrode  is  assumed  to  be  spherical  particles  of  electron-  and 
ion-conducting  materials  packed  together  at  random.  Due  to  ran¬ 
domness  of  these  models,  the  random  resistor  network  models 
and  random  packing  sphere  models  are  also  referred  to  as  Monte 
Carlo  resistor  network  models  and  Monte  Carlo  packing  sphere 
models,  respectively.  On  the  other  hand,  macroscopic  porous- 
electrode  models  are  based  on  the  assumption  that  the  composite 
electrode  is  represented  by  the  particles  of  electron-  and  ion¬ 
conducting  materials  packed  together  at  random  and  disregards 
the  actual  geometric  details  of  the  individual  particles.  Instead,  the 


composite  electrode  is  described  in  terms  of  continuous-average 
quantities.  Hence,  in  the  present  model,  the  macroscopic  porous- 
electrode  approach  is  used  to  model  reaction  zone  layers,  which 
according  to  Sunde  [15]  is  the  most  practical  method  for  applied 
research. 

The  processes  that  need  to  be  modeled  in  the  reaction  zone  lay¬ 
ers  are  transport  of  multi-component  mixture  of  species  to  and 
from  the  reaction  sites  along  with  chemical  (water-gas  shift  and 
methane  reforming)  and  electrochemical  reactions,  respectively 
transport  of  electrons  and  ions  in  the  electron-  and  ion-conducting 
particles,  and  transfer  of  energy  due  to  heat  conduction,  species 
diffusion,  heat  generation/consumption  due  to  chemical  and  elec¬ 
trochemical  reactions  and  electron  and  ion  migration. 

The  conservation  of  species  equation  governing  the  transport 
of  multi-component  mixture  in  the  reaction  zone  layers  is  similar 
to  the  backing  layers  given  in  Eq.  (1)  except  the  source  term, 
which  accounts  for  the  consumption  or  production  of  species 
due  to  electrochemical  reactions  in  addition  to  the  chemical 
reactions  described  in  the  anode  backing  layer.  Again,  the  modified 
Stefan-Maxwell  equations  given  in  Eq.  (4)  are  used  to  model  the 
diffusion  flux  in  the  reaction  zone  layers.  The  equations  governing 
the  transport  of  electrons  and  oxide  ions  in  the  electron-  and 
ion-conducting  particles  of  the  reaction  zone  layers  are  described 


as  follows; 

Anode  reaction  zone  layer: 

Electronic  charge  :  V  •  ~fe  =  Jla  (30) 

Ionic  charge  :  V  •  ~f  {  =  -JZa  (31) 

Cathode  reaction  zone  layer: 

Electronic  charge  :  V  •  ~fe  =  -Kc  (32) 

Ionic  charge:  V  •  =  7£c  (33) 


where  JZa  and  JZC  are  the  volumetric  current  densities  produced  in 
the  anode  and  cathode  reaction  zone  layers  due  to  H2  oxidation  and 
02  reduction  reactions,  respectively,  which  are  represented  by  the 
general  Butler-Volmer  equations,  expressed  as 


Ka=A„jfr 


ref 


Ch7 


Ch- 


2,  re/ 


Kc=AvJ°>ef 


COy 


Co- 


'2, ref 


Y\\y 


YOn 


(34) 


h(-?)-'(T)) 


(35) 


where  Av  is  the  actual  reactive  surface  area  per  unit  volume, 
and  J^rejr  are  the  reference  exchange  current  densities  for  H2  oxida¬ 
tion  and  02  reduction  reactions  at  reference  concentrations  Ch2  ref 
and  cq2  re/,  respectively,  a  is  the  charge  transfer  coefficient  (or  sym¬ 
metry  factor),  whose  value  lies  between  zero  and  one,  n  is  the 
number  of  electrons  participating  in  the  electrochemical  reaction 
and  ri  is  the  electrode  activation  overpotential  (anode  or  cathode), 
defined  as 


r)  =  <pi-  4>e 


(36) 
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In  order  to  enhance  the  predictive  capability  of  the  model,  the 
expression  used  to  model  the  reactive  surface  area  per  unit  volume 
is  based  on  random  packing  of  binary  spherical  particles  developed 
by  Costamagna  and  Honegger  [41  ],  which  is  given  as 


Av  =  jrsin2  0r2elntnelnio^^-pelpio 


(37) 


where  0  is  the  contact  angle  between  the  electron-  and  ion¬ 
conducting  particles  in  the  reaction  zone  layers,  rd  is  the  radius  of 
the  electron-conducting  particles,  nt  is  the  total  number  of  particles 
per  unit  volume,  nd  and  nio  are  the  number  fractions  of  the  electron- 
and  ion-conducting  particles  in  the  reaction  zone  layers,  respec¬ 
tively,  Zd  and  Zio  are  the  coordination  numbers  of  the  electron-  and 
ion-conducting  particles  in  the  reaction  zone  layers,  respectively, 
Z  is  the  total  average  number  of  contacts  of  each  particle  and  pd 
and  pio  are  the  probabilities  of  the  electron-  and  ion-conducting 
particles  in  the  reaction  zone  layers,  respectively. 

The  parameters  required  to  obtain  the  reactive  surface  area  per 
unit  volume  (Av)  are  calculated  using  the  following  expressions 
[41,15]: 

nt  = - - =-  (38) 

(4/3 )jt r2  [ne/  +  (1  -  ne,)(rio/rely] 

where  e  is  the  porosity  of  the  reaction  zone  layers: 

0 

nel  = 


(39) 


[<£>  +  ((! -<2>)/(rioAe()3)] 
where  0  is  the  volume  fraction  of  the  electron-conducting  particles 
in  the  reaction  zone  layers: 

(Z  —  3) 


Ze\  —  3  + 


Zio  —  3  + 


[^e/T(l  ^el)iXio/^el) 
(Z  -  3  )(ri0/re,)2 


(40) 


(41) 


[ftel  T  (1  ^e/XHo/^e/)  ] 

where  Z  is  the  total  average  coordination  number,  equal  to  6  [15]: 

(42) 

Pio  =  [1  —  (2  -  Z,-0_;0/2)2'5]0  4  (43) 

where 


Pet  =  [1  —  (2  —  Ze|_e;/2)2'5]0  4 

0.4 


Zel-el  ~ 


ne\Z2d 

z 

riinZl 


(44) 


(45) 


The  species  source  terms  in  the  anode  and  cathode  reaction  zone 
layers  due  to  electrochemical  reactions  are  related  to  volumetric 
current  density  produced  through  Faraday’s  law  of  electrochemical 
reaction,  and  are  described  as  follows: 

Anode  reaction  zone  layer: 


b  2^a 

Ss^  =  —^i ~ 


(46) 


where  7 Za  is  the  volumetric  current  density  produced  due  to  H2 
oxidation  reaction  given  in  Eq.  (34),  n  is  the  number  of  electrons 
participating  in  the  electrochemical  reaction,  Fis  the  Faraday’s  con¬ 
stant,  and  uh2  is  the  stoichiometric  coefficient  of  the  H2  involved  in 
the  oxidation  reaction  expressed  in  the  following  form  [26]: 

H2  +  02“  -  H20  2e~  (47) 


Cathode  reaction  zone  layer: 


Ss,o2  =  - 


^o2  'Kc 

nF 


(48) 


where  7 Zc  is  the  volumetric  current  density  produced  due  to  02 
reduction  reaction  given  in  Eq.  (35),  and  vq2  is  the  stoichiometric 


coefficient  of  the  02  involved  in  the  reduction  reaction  expressed  in 
the  following  form  [26]: 

lo2-02-^2e-  (49) 

Expressing  the  electronic  and  ionic  current  densities  in  terms  of 
electronic  and  ionic  potentials  through  Ohm’s  law,  Eqs.  (30)-(33) 


become 

Anode  reaction  zone  layer: 

Electronic  charge  :  V  •  (a^ffV0e)  =  -7 Za  (50) 

Ionic  charge  :  V  •  =  7 Za  (51) 

Cathode  reaction  zone  layer: 

Electronic  charge  :  V  •  (a^ffV0e)  =  7£c  (52) 

Ionic  charge  :  V  •  (/c^ffV0j)  =  -7 Zc  (53) 

where  the  effective  electronic  and  ionic  conductivities  in  the  reac¬ 
tion  zone  layers  are  determined  as  follows: 

a*=*(^)a  (54) 

icf  =  {  (55) 


where  0  is  the  volume  fraction  of  the  electron-conducting  particles 
in  the  reaction  zone  layers,  a  and  /c  are  the  conductivities  of  pure 
electron-  and  ion-conducting  materials,  respectively. 

The  conservation  of  energy  equation  in  the  reaction  zone  layers 
is  similar  to  the  backing  layers  except  the  energy  source  term,  which 
is  expressed  as 


V  •  (kfVT) 


■VT- 


+  Se  —  0 


(56) 


where  fc®ff  is  the  effective  thermal  conductivity  of  the  reaction  zone 
layer,  calculated  as 

kf  =  skj  +  (1  -  s)[0kel  +  (1  -  0)kio]  (57) 

where  kf  is  the  thermal  conductivity  of  the  fluid  mixture  in  the 
reaction  zone  layers  obtained  using  Eq.  (26),  0  is  the  volume  frac¬ 
tion  of  the  electron-conducting  particles  in  the  reaction  zone  layers, 
kd  and  ki0  are  the  thermal  conductivities  of  pure  electron-  and 
ion-conducting  materials,  respectively. 

The  energy  source  term  Se  in  the  reaction  zone  layers  is  repre¬ 
sented  as 

Anode  reaction  zone  layer: 


eff  ^  ^eff 
rl 


7 

nF 


(— r(As)H2,0X)  +  Ka/? 


(58) 


where  the  first  term  on  the  right-hand  side  (RHS),  ([/2l)/(o’Jjff), 
represents  the  Joule  heating  effect  due  to  the  resistance  to  the 
electron  flow  in  the  electron-conducting  particles  of  the  anode 
reaction  zone  layer,  the  second  term,  ( L/?  I  )/(^ff),  represents  the 
Joule  heating  effect  due  to  the  resistance  to  the  ion  flow  in  the 
ion-conducting  particles  of  the  anode  reaction  zone  layer,  the  third 
term,  (|7£fl/nF|(-r(As)H  ox)),  represents  the  reversible  heat  gen¬ 
eration  due  to  the  electrochemical  H2  oxidation  reaction,  and  the 
last  term,  7 Zar),  represents  the  irreversible  heat  generation  due  to 
electrochemical  H2  oxidation  reaction. 

In  order  to  find  the  reversible  heat  generation  due  to  the  elec¬ 
trochemical  reaction  in  the  anode  reaction  zone  layer,  the  entropy 
change  between  the  reactants  and  the  products  has  to  be  deter¬ 
mined.  However,  the  entropy  of  formation  for  the  ions  and  electrons 
are  not  known.  Therefore,  the  entropy  change  for  an  overall  SOFC 
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Table  2 

Governing  equations  in  the  reaction  zone  layers. 


Governing  equations 

Anode  reaction  zone  layer 

Cathode  reaction  zone  layer 

Species 

v-7r,=ss>( 

V-7?,=Ss,o2 

Modified  S-M 

=  Ej.i 

Electronic  Charge 

V  ■  (a?V<j>e)  =  -Ra 

V  •  (<T'ffV0e)  =  ftc 

Ionic  Charge 

V  •  V*)  =  TZa 

V  ■  V0r)  =  -Tic 

Energy 

v'(kTvf-(El^)'VT+^  +  §  + 

v-(*fvr)-(X)^c„H;)-vr  +  ^  + 

|Sj  (-T?h2(As)„2)  +  K.,- (X^h.Sv)  =0 

1  w  1  (—7(1  ~  ?h2  X  As)H2  )  +  TZcH  =  0 

reaction  is  determined  in  place  of  overall  half-cell  reactions.  The 
overall  SOFC  reaction  is 

H2  +  \o2  =  H20(gj  (59) 

The  change  in  entropy  between  the  products  and  reactants  for 
the  above  overall  SOFC  reaction  is 

(  As)h2  =  5h2o  -  Sh2  -  2^o2  (60) 

The  entropy  change  in  the  anode  reaction  zone  layer  during  the 
H2  oxidation  reaction  is  assumed  to  be  a  fraction  of  the  entropy 
change  for  the  overall  SOFC  reaction,  so  that 


Se  = 


j/?i ,  j/?i 

rjeff  ^  „eff 
°rl  Krl 


JZa 

nF 


(-7?h2(As)H2)  +  JZarj 


(61) 


where  §h2  is  the  fraction  of  the  entropy  change  for  the  overall  SOFC 
reaction  of  H2. 

Cathode  reaction  zone  layer: 


Se 


_[ £1,1/?!, 

n-eff  ^  ^eff  _1" 
°rl  Krl 


nF 


[ — T(1  -  §h2)(As)h2]  +  7lcri 


(62) 


where  the  first  term  on  the  right-hand  side,  (|/|  l)/(^1ff),  repre¬ 
sents  the  Joule  heating  effect  due  to  the  resistance  to  the  electron 
flow  in  the  electron-conducting  particles  of  the  cathode  reac¬ 
tion  zone  layer,  the  second  term,  (\jf\)/{icff),  represents  the  Joule 
heating  effect  due  to  the  resistance  to  the  ion  flow  in  the  ion¬ 
conducting  particles  of  the  cathode  reaction  zone  layer,  the  third 
term,  {\Kc/nF\(-T(l  -  £h2)(As)h2 >ox)),  represents  the  reversible 
heat  generation  due  to  the  electrochemical  02  reduction  reaction, 
and  the  last  term,  Kcr),  represents  the  irreversible  heat  generation 
due  to  the  electrochemical  02  reduction  reaction.  The  governing 
equations  in  the  reaction  zone  layers  are  summarized  in  Table  2. 


2.3.  Electrolyte  layer 

The  electrolyte  layer  in  SOFC  is  fully  dense  with  no  intercon¬ 
nected  porosity.  The  function  of  the  electrolyte  layer  is  to  conduct 
oxide  ions  produced  in  the  cathode  reaction  zone  layer  to  the  anode 
reaction  zone  layer,  thus  completing  the  electrical  circuit.  The  pro¬ 
cesses  which  need  to  be  modeled  are  the  transport  of  oxide  ions 
and  transfer  of  energy  due  to  heat  conduction  and  heat  genera¬ 
tion  due  to  ion  migration.  The  conservation  equation  governing  the 
processes  in  the  electrolyte  layer  are  given  as 

Ionic  charge  :  V  •  ~f  {  =  0  (63) 

Energy  :  V  if  =  Se  (64) 

where  is  the  ionic  current  density,  which  is  equal  to  the  cell 
current  density  and  Se  is  the  energy  source  term  due  to  the 
resistance  to  the  transport  of  oxide  ions,  which  is  expressed  as 

I2 

Se  =  J-  (65) 

where  k  is  the  ionic  conductivity  of  the  electrolyte. 


Using  Ohm’s  law  and  Fourier’s  law  of  heat  conduction,  the  ionic 
current  density  and  heat  flux  are  expressed  in  terms  of  the  ionic 
potential  and  temperature,  respectively;  the  above  conservation 
equations  become 

Ionic  charge  :  V  •  (/cV0 f)  =  0  (66) 

Energy:  V  •  (kdeVT)  +  J1  =  0  (67) 

where  k  and  kde  are  the  ionic  and  thermal  conductivities  of  the 
electrolyte,  respectively.  The  governing  equations  in  the  electrolyte 
layer  are  summarized  in  Table  3. 

2.4.  Boundary  conditions 

Since  SOFC  modeling  involves  interdependent  transport  pro¬ 
cesses  in  different  components,  the  specification  of  boundary 
conditions  includes  both  external  and  internal  boundary  condi¬ 
tions  [27].  The  locations  at  which  boundary  conditions  are  required 
in  the  model  are  illustrated  in  the  schematic  shown  in  Fig.  1.  The 
locations  of  the  external  boundary  conditions  are  represented  by 
E1-E6;  whereas,  the  locations  of  the  internal  boundary  conditions 
are  represented  by  11-14.  Since  a  planar  SOFC  is  symmetric  about 
the  middle  of  the  flow  channel,  symmetric  boundary  conditions 
are  applied  at  El.  The  boundaries  E2  and  E6  represent  the  inter¬ 
faces  between  the  gas  flow  channels  and  the  backing  layers,  where 
composition  of  gaseous  species  are  specified,  electronic  current 
densities  are  insulated,  and  the  temperature  at  these  boundaries 
is  specified  as  the  inlet  temperatures  of  gas  streams.  The  bound¬ 
aries  E3  and  E5  represent  the  land  portions  of  the  interconnects 
which  are  in  contact  with  the  backing  layers,  where  the  mass  fluxes 
of  species  are  zero,  electronic  potential  is  set  as  ground  potential  at 
boundary  E3,  and  either  cell  potential  or  current  density  is  specified 
at  boundary  E5  depending  upon  the  approach  used  in  the  particu¬ 
lar  simulation.  Since  the  computational  domain  only  includes  the 
land  portions  of  the  interconnects,  the  temperatures  are  specified 
at  boundaries  E3  and  E5  as  first  approximations,  instead  of  con¬ 
tinuous  heat  flux  boundary  conditions,  which  are  applicable  when 
the  interconnects  are  included  in  the  computational  domain.  The 
boundary  E4  represents  the  external  boundary,  where  insulated 
boundary  conditions  are  applicable,  which  means  the  molar  flux 
of  species,  current  densities  (electronic  and  ionic),  and  heat  flux 
are  zero. 

The  internal  boundaries  II  and  14  represent  the  interfaces 
between  the  backing  and  reaction  zone  layers.  The  boundary  condi¬ 
tions  at  these  boundaries  are  continuous  flux  boundary  conditions 
for  the  species  and  electronic  potential,  where  the  molar  flux 

Table  3 

Governing  equations  in  the  electrolyte  layer. 
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Table  4 

Boundary  conditions. 


Boundary/internal  locations  Conditions 

El 

3x,- 

_  d<pe  _  _  dT 

=  0 

dy 

~  dy  ~  3 y~  dy 

E2  and  E6 

Xi  : 

=  specified;  je  ■  it 

=  0;  T  =  specified 

E3 

aT 

■  It  =  0;(pe  =  0;T 

=  specified 

E4 

AT 

■ft  =  Te- It  =7 

r  -ft  =  t  ■  Tt  = 

0 

E5 

aT 

■  Tt  =  0;  7e  ■  ft  = 

J'.T  =  specified 

11  and  14 

aT 

•  ft'lw  =  'tff  •  Ttlri; 

7  e  '  1  bl  =  7  e 

Tflrl 

r 

•Tf|  =  0;t  •  ■JflM 

=  t  Tt\rt 

12  and  13 

aT 

•  Tt  =  0;  Te  •  Tt  = 

0-,7fTt=J 

t 

■Tt\rl  =  t  Tt\el 

and  electronic  current  density  are  continuous;  whereas,  insulated 
boundary  condition  applies  for  the  ionic  potential,  which  implies 
ionic  current  density  is  zero.  The  interfaces  between  the  reaction 
zone  layers  and  the  electrolyte  layer  are  represented  by  12  and 
13.  Since  the  electrolyte  is  an  electron  insulator  with  no  intercon¬ 
nected  porosity,  the  molar  flux  of  species  and  electronic  current 
density  are  zero  at  these  interface  boundaries.  Whereas,  the  bound¬ 
ary  condition  for  the  ionic  potential  depends  on  the  approach  used 
in  fuel  cell  modeling.  As  it  is  known  in  fuel  cell  modeling,  either 
current  density  is  specified  to  obtain  the  cell  potential  or  vice 
versa.  When  current  density  is  specified,  the  ionic  current  den¬ 
sity  at  these  boundaries  (12  and  13)  is  set  equal  to  the  cell  current 
density;  however,  when  cell  potential  is  specified  to  obtain  the  cur¬ 
rent  density,  the  ionic  current  density  at  these  boundaries  (12  and 
13)  is  continuous.  Similarly,  the  heat  flux  is  continuous  at  these 
boundaries  irrespective  of  the  approaches  used.  Finally,  the  bound¬ 
ary  conditions  expressed  in  mathematical  form  are  summarized  in 
Table  4. 

3.  Numerical  implementation 

The  governing  equations  tabulated  in  Tables  1-3  need  to  be 
discretized  to  obtain  the  numerical  solution.  Since  the  mathemat¬ 
ical  model  governing  the  processes  in  different  layers  of  an  SOFC 
is  based  on  the  conservation  principle,  the  discretization  scheme 
commonly  used  for  problems  involving  conservative  laws  is  a  well- 
known  method,  often  referred  to  as  finite  volume  method  or  control 
volume  method  [28,29].  The  application  of  finite  volume  method  to 
fuel  cell  related  problems  has  already  been  recognized  by  various 
researchers  [30,31]. 

The  governing  discretized  equations  form  a  system  of  algebraic 
equations,  whose  solution  can  be  obtained  using  a  direct  or  an  iter¬ 
ative  method.  For  the  present  problem,  which  is  highly  non-linear, 
involving  interdependent  variables  with  multi-component  mixture 
of  species,  an  iterative  solver  is  used  to  obtain  the  solution  of  the 
discretized  governing  equations. 

The  grid  (numerical  mesh)  is  uniform  locally  but  non-uniform 
globally,  which  means  the  mesh  size  is  uniform  within  a  layer  or 
component  of  the  computational  domain.  Finer  mesh  sizes  are  used 
in  the  reaction  zone  layers  than  in  the  backing  layers  to  capture 
the  electrochemical  reactions  and  other  processes.  Grid  sensitiv¬ 
ity  tests  are  performed  to  ensure  grid  independent  solution.  For 
instance,  beyond  the  grid  size  of  150  x  50,  no  significant  change 
is  observed  in  the  species  distribution  in  different  layers  of  the 
cell. 

4.  Model  validation 

Model  validation  is  as  important  as  model  development.  It  helps 
in  determining  the  range  of  validity  and  accuracy  of  the  model. 
However,  only  few  researchers  [4,5,32]  have  shown  the  validation 
of  their  SOFC  models.  The  parameters  used  in  the  validation  of  the 
present  model  are  listed  in  Table  5.  Most  of  the  parameters  are 


Table  5 

Parameters  used  for  model  validation. 


Source :  Yakabe  et  al.  [4]. 

obtained  from  Yakabe  et  al.  [4],  who  reported  the  measured  data  for 
concentration  overpotential  for  a  ternary  mixture  at  different  cur¬ 
rent  densities.  The  thickness  of  the  anode  is  reported  to  be  2  mm,  of 
which  50  |jim  is  treated  as  the  thickness  of  the  anode  reaction  zone 
layer,  which  is  in  accordance  with  the  depth  of  the  reaction  zone 
reported  by  Lehnert  et  al.  [11  ]  for  an  anode  thickness  of  2  mm. 

Fig.  2  shows  the  comparison  between  the  model  prediction  and 
measured  data  at  0.7  A  cm-2.  The  abscissa  of  Fig.  2  represents  the 
variation  in  reactant  concentration,  while  the  ordinate  represents 
the  difference  between  the  actual  and  reference  concentration  over¬ 
potential.  The  reference  concentration  overpotential  was  measured 
when  there  was  no  argon  in  the  system  at  hydrogen  mole  fraction 
equal  to  0.8.  Argon  was  added  to  the  system  to  vary  the  reactant  con¬ 
centration  such  that  the  ratio  of  the  mole  fractions  of  hydrogen  and 
water  vapor  is  4:1.  It  can  be  seen  that  increasing  the  reactant  con¬ 
centration  decreases  the  concentration  overpotential  difference. 


Fig.  2.  Model  validation. 
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Total  pressure,  p  (atm) 

1.0 
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Anode  electrode  layer  thickness,  ta  (pm) 

1950.0 
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(pan) 

Resistivity  of  electron-conducting 

2.98  x  10-5  exp  ( - )  [33] 
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Resistivity  of  ion-conducting  paricles,  pi0 
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Porosity,  e 
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Tortuosity,  r 

4.5 
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15°  [41] 
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Radius  of  electron-conducting  particles,  re/ 
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Radius  of  ion-conducting  particles,  rio  (pm) 

0.1  [13,15] 

Volume  fraction  of  electron-conducting 

0.5  [12,15] 

particles,  0 

Reference  H2  concentration,  cH2 

10.78 

(mole  m-3) 

Reference  exchange  current  density  for  H2 
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oxidation  Jjh  (Am-2) 

Reaction  order  for  H2  oxidation,  yH2 
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Table  6 

Fuel  composition  used  in  the  simulation. 


Source :  Lehnert  et  al.  [11  ]. 

This  is  due  to  increased  reactant  concentration  at  the  reactant  sites, 
which  results  in  reducing  the  concentration  overpotential.  Further, 
it  can  be  seen  that  the  concentration  overpotential  difference  pre¬ 
dicted  by  the  model  is  in  excellent  agreement  with  measured  data. 
The  coefficient  of  determination  or  R2  value,  an  indicator  between 
0  and  1,  reveals  how  closely  the  predicted  values  correspond  to 
measured  data,  is  obtained  as  0.989. 

5.  Results  and  discussion 

The  developed  model  not  only  predicts  the  concentration  over- 
potential  but  also  predicts  the  activation  and  ohmic  overpotentials, 
which  in  turn  helps  in  predicting  the  overall  performance  of  an 
anode-supported  SOFC.  The  fuel  composition  and  the  base  case 
parameters  used  in  the  simulation  are  listed  in  Tables  6  and  7, 
respectively. 

Fig.  3  shows  the  cell  performance  resulting  from  the  fuel  compo¬ 
sition  and  base  case  parameters  listed  in  Tables  6  and  7,  respectively. 
The  abscissa  of  Fig.  3  represents  the  current  density  (load),  whereas 
the  primary  (left)  and  secondary  (right)  ordinate  represent  the 
cell  potential  and  power  density,  respectively.  The  solid  line  cor¬ 
responds  to  the  cell  potential  and  power  density  when  all  the 
overpotentials  are  subtracted  from  the  reversible  cell  potential.  The 
overpotentials  include  anode  overpotential,  cathode  overpotential 
and  electrolyte  overpotential.  Anode  and  cathode  overpotentials 
include  activation  overpotentials  due  to  the  resistance  to  the  charge 
transfer  reactions,  ohmic  overpotentials  due  to  the  resistance  to 
the  flow  of  electrons  and  ions  in  the  reaction  zone  layers  and 
the  resistance  to  the  flow  of  electrons  in  the  backing  layers,  and 
concentration  overpotentials  due  to  the  resistance  to  the  flow  of 
reactant  species  through  the  void  spaces.  Conversely,  the  dashed 
lines  correspond  to  the  cell  potential  when  anode  overpoten¬ 
tial,  electrolyte  overpotential,  and  cathode  overpotential  are  not 


Fig.  3.  Base  case  performance  of  anode-supported  an  SOFC. 


Table  7 

Base  case  parameters  used  in  the  simulation. 


included  in  obtaining  the  polarization  curves.  It  can  be  seen  that 
the  solid  line  representing  the  actual  overpotential  exhibits  regions 
of  activation  and  concentration  overpotentials  at  low  and  high  cur¬ 
rent  densities,  respectively.  Specifically,  these  regions  of  activation 
and  concentration  overpotentials  are  due  to  the  cathode  activation 
and  anode  concentration  overpotentials,  which  is  evident  from  the 
dashed  lines  neglecting  cathode  and  anode  overpotentials,  respec¬ 
tively.  Additionally,  it  is  observed  that  the  maximum  power  density 
obtained  for  an  anode-supported  SOFC  using  the  fuel  composi¬ 
tion  and  base  case  parameters  listed  in  Tables  6  and  7  is  about 
0.4  Wcrrr2.  It  is  worthwhile  to  note  that  the  average  values  of 
activation  overpotentials  are  used  in  obtaining  the  polarization 
curves,  since  the  values  of  activation  overpotentials  vary  within  the 
reaction  zone  layers  with  maximum  at  the  interface  between  the 
reaction  zone  and  electrolyte  layers  and  minimum  at  the  interface 
between  the  backing  and  reaction  zone  layers. 

It  is  also  seen  from  Fig.  3  that  the  profiles  without  the  anode 
and  cathode  overpotentials  exhibit  almost  similar  performances 
for  current  density  upto  about  lAcm-2,  indicating  comparable 
magnitudes  of  anode  and  cathode  overpotentials.  Flowever,  at  cur¬ 
rent  densities  beyond  lAcrrr2,  the  profile  without  the  anode 
overpotential  shows  better  performance.  The  reason  for  better 
performance  at  higher  current  densities  is  because  of  anode  con¬ 
centration  overpotential,  which  is  not  included  in  the  profile 
without  anode  overpotential.  Further,  it  is  observed  that  the  con¬ 
tribution  of  electrolyte  overpotential  is  not  negligible  even  though 
its  thickness  is  negligible  relative  to  the  thickness  of  the  anode.  To 
better  understand  these  observations,  the  anode  overpotentials  at 
the  base  case  condition  are  further  investigated,  and  a  typical  result 
is  shown  in  Fig.  4.  It  is  clear  that  the  ohmic  overpotential  is  the  sin¬ 
gle  largest  contributor  to  the  overall  anode  overpotential  and  the 
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Fig.  4.  The  breakdown  of  the  anode  overpotentials  at  the  base  case  conditions  given 
in  Table  7. 


contribution  of  the  activation  overpotential  is  negligible.  Moreover, 
the  concentration  overpotential  in  the  thick  anode  remains  almost 
constant  and  close  to  zero  for  most  of  the  current  density  range 
considered  in  the  present  simulation  before  starting  to  increase 
beyond  about  1.1  A/cm2.  The  reasons  for  the  anode  concentration 
overpotential  being  negligible  even  in  the  thick  anode  are  because 
of  low  fuel  utilization  (10%)  and  high  reactant  concentration  due  to 
the  chemical  reactions  (reforming  and  water-gas  shift  reactions). 
For  instance,  methane  reforming  reaction  produces  three  moles 
of  H2for  every  mole  of  CFUconsumed  during  the  reaction;  simi¬ 
larly,  water-gas  shift  reaction  produces  a  mole  of  H2  for  every  mole 
of  CO  consumed  during  the  reaction.  A  useful  comparison  of  the 
anode  concentration  overpotential  with  and  without  considering 
these  chemical  reactions  is  shown  in  Fig.  5.  Even  at  low  fuel  utiliza¬ 
tion,  the  concentration  overpotential  in  the  anode-supported  SOFC 
becomes  significant  when  there  are  no  chemical  reactions  in  the 


anode  and  can  be  as  high  as  0.1  V  at  higher  current  densities.  How- 
ever,  it  is  still  three  orders  of  magnitude  smaller  than  the  anode 
ohmic  overpotential. 

In  contrast  to  the  anode  concentration  overpotential,  the  cath¬ 
ode  concentration  overpotential  in  an  anode-supported  SOFC  is 
negligible.  This  is  again  due  to  the  small  thickness  of  the  cathode 
and  low  oxidant  utilization  (10%).  However,  the  contributions  of 
the  cathode  activation  and  cathode  ohmic  overpotentials  are  sig¬ 
nificant  to  the  total  cell  potential  loss,  which  is  evident  from  Fig.  6. 
Further,  it  is  clear  from  Figs.  4  and  6  that  the  ohmic  overpotential  is 
the  single  largest  contributor  to  the  cell  potential  loss  even  in  the 
anode-supported  SOFCs,  which  is  often  reported  to  be  negligible. 
Furthermore,  it  can  be  seen  from  Figs.  4  and  6  that  the  profiles  of 
the  anode  and  cathode  ohmic  overpotentials  show  different  degrees 
of  non-linear  trends  with  current  density.  This  is  due  to  the  larger 
thickness  of  the  anode  reaction  zone  layer  than  the  cathode  reaction 
zone  layer  resulting  in  the  spatial  limitation  of  oxide  ion  trans¬ 
port  in  the  anode  reaction  zone  layer  [40].  Further,  the  temperature 
dependent  ionic  conductivity  of  the  ion-conducting  particles  in  the 
reaction  zone  layers  contributes  to  the  non-linear  behavior  of  the 
ohmic  overpotentials. 

In  order  to  probe  the  robustness  of  the  model,  a  parametric 
study  has  been  performed  to  examine  the  effect  of  operating  and 
design  conditions  on  the  performance  of  an  anode-supported  SOFC. 
Moreover,  as  it  is  rightly  stated  by  Alkhateeb  et  al.  [35]  that  a  com¬ 
prehensive  model  would  analyze  the  effect  of  macroscopic  and 
microscopic  characteristics  of  the  electrodes  on  the  performance 
of  the  cell.  To  start  with,  the  effect  of  operating  temperature  on  the 
performance  of  an  anode-supported  SOFC  is  shown  in  Fig.  7.  Anode- 
supported  SOFCs  typically  operate  between  823  I<  and  1073  K[3].In 
the  present  simulation,  temperature  in  the  typical  range  is  varied 
to  examine  its  effect  on  the  performance  of  an  anode-supported 
SOFC.  The  reaction  rate  constants  for  the  chemical  reactions  are 
valid  for  the  temperature  range  considered  in  the  present  sim¬ 
ulation  [20].  Further,  the  other  operating  and  design  parameters 
are  kept  constant  in  accordance  with  the  base  case  parameters.  It 
can  be  seen  that  increasing  the  operating  temperature  of  the  cell 
increases  the  performance  of  an  anode-supported  SOFC.  It  can  also 
be  seen  that  increasing  the  operating  temperature  increases  the 
limiting  current  density,  which  corresponds  to  the  current  density 
at  the  zero  cell  potential.  Moreover,  the  power  density  increases 


Current  Density  (A/cm2) 

Fig.  5.  Comparison  of  the  anode  concentration  overpotential  at  the  base  case  con¬ 
ditions  given  in  Table  7.  The  chemical  reactions  refer  to  the  methane  reforming  and 
water-gas  shift  reactions  that  produce  hydrogen  in  the  anode. 


0.2 


>  0.15 

.2 

c 

a> 

o 

a 

S  o.i 

> 

o 

<D 

•a 

o 

.c 

Q  0.05 


0 


0  0.2  0.4  0.6  0.8  1  1.2  1.4  1.6 

Current  Density  (A/cm2) 

Fig.  6.  The  breakdown  of  theCathode  overpotentials  at  the  base  case  conditions 
given  in  Table  7. 
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Fig.  7.  Effect  of  operating  temperature  on  the  performance  of  an  anode-supported 
SOFC. 

with  the  increase  in  operating  temperature,  with  its  peak  shifted 
towards  higher  current  densities  at  higher  temperatures.  The  rea¬ 
son  for  the  increase  in  the  cell  performance  with  temperature  is 
primarily  due  to  the  increase  in  the  temperature  dependent  ionic 
conductivity  of  the  ion-conducting  particles  in  the  reaction  zone 
layers,  which  in  turn  reduces  the  contribution  of  the  ohmic  over¬ 
potential  of  the  electrodes.  Further,  the  ionic  conductivity  of  the 
electrolyte  increases  with  the  increase  in  operating  temperature, 
resulting  in  better  cell  performance.  The  increase  in  the  limiting 
current  density  with  operating  temperature  is  due  to  the  fact  that 
the  molecular  diffusivity  of  the  species  increases  with  tempera¬ 
ture,  which  reduces  the  resistance  to  the  mass  transport  in  the  thick 
anode,  and  thereby  reduces  the  anode  concentration  overpotential 
with  temperature.  Although  there  are  many  incentives  in  reducing 
the  operating  temperature  but  there  is  a  significant  drop  in  cell  per¬ 
formance  just  by  reducing  the  temperature  from  1073  K  to  973  K, 
which  is  evident  from  Fig.  7.  Flence,  the  ionic  conductivity  of  the 
ion-conducting  particles  in  the  reaction  zone  layers  and  electrolyte 
need  to  be  enhanced  in  order  to  operate  anode-supported  SOFCs 
below  1073  K. 

The  thickness  of  the  reaction  zone  layers  is  one  of  the  important 
parameters  in  SOFC  electrodes.  Different  values  are  reported  in  lit¬ 
erature,  based  on  the  thickness  of  the  electrode  ranging  from  10  [xm 
to  50  pan  [11,17,36-39].  The  effect  of  the  anode  reaction  zone  thick¬ 
ness  on  the  performance  of  an  anode-supported  SOFC  is  shown  in 
Fig.  8.  The  operating  and  other  design  parameters  are  kept  the  same 
as  those  in  Table  7.  In  addition,  the  combined  thickness  of  the  anode 
backing  layer  and  the  reaction  zone  layer  is  kept  at  2  mm.  Increasing 
the  reaction  zone  thickness  increases  TPBs  in  the  reaction  zone  layer 
resulting  in  the  increased  rate  of  electrochemical  reaction,  which 
in  turn  reduces  the  activation  overpotential.  However,  increasing 
the  reaction  zone  thickness  increases  the  distance  through  which 
oxide  ions  and  electrons  migrate  to  reach  the  reaction  sites.  Due  to 
the  poor  ionic  conductivity  of  the  ion-conducting  particles,  increas¬ 
ing  the  reaction  zone  thickness  increases  the  ohmic  overpotential. 
The  balance  between  the  decrease  in  the  activation  overpotential 
and  the  increase  in  the  ohmic  overpotential  reflects  the  improve¬ 
ment  in  the  cell  performance  with  respect  to  the  increase  in  the 
anode  reaction  zone  thickness.  A  close  observation  of  Fig.  8  reveals 
that  there  is  no  significant  increase  in  the  performance  when  the 
anode  reaction  zone  thickness  is  increased  from  30  |xm  to  70  |xm, 
in  fact,  the  performance  is  reduced  at  higher  current  densities. 

The  effect  of  porosity  on  the  performance  of  an  anode-supported 
SOFC  is  illustrated  in  Fig.  9.  All  other  operating  and  design  param¬ 


Fig.  8.  Effect  of  the  anode  reaction  zone  thickness  on  the  performance  of  an  anode- 
supported  SOFC. 

eters  are  kept  the  same  as  for  the  condition  shown  in  Table  7.  It 
is  seen  that  increasing  the  porosity  of  the  porous  layers  decreases 
the  performance  of  an  anode-supported  SOFC;  however,  the  per¬ 
formance  increases  with  increasing  the  porosity  at  higher  current 
densities.  This  is  due  to  the  reduction  in  mass  transport  resistance 
with  increasing  porosity  of  the  porous  layers  at  higher  current  den¬ 
sities.  In  contrary,  increasing  the  porosity  decreases  the  effective 
conductivities  of  the  porous  layers  resulting  in  the  increased  con¬ 
tribution  of  ohmic  overpotentials  for  most  of  the  current  density 
range  considered  in  the  present  simulation. 

Fig.  10  shows  the  effect  of  tortuosity  of  the  porous  layers  on 
the  performance  of  an  anode-supported  SOFC.  Again,  all  the  other 
design  and  operating  parameters  are  set  equal  to  the  base  case 
parameters  given  in  Table  7.  One  of  the  primary  deficiencies  of  some 
of  the  earlier  models  is  the  need  to  invoke  unreasonably  large  tortu¬ 
osity  values  in  the  range  of  10-17  [5,41  ].  The  reason  of  invoking  such 
large  tortuosities  is  to  produce  the  concentration  overpotential  [40]. 
However,  tortuosities  for  the  porous  layers  of  SOFC  have  been  exper¬ 
imentally  determined  to  be  in  the  range  of  2.5-6.0  [11,17,40,42]. 
Therefore,  in  the  present  simulation,  tortuosity  of  the  porous  layers 
has  been  varied  between  3.0  and  6.0.  Increasing  the  tortuosity  of 


Fig.  9.  Effect  of  porosity  on  the  performance  of  an  anode-supported  SOFC. 
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Fig.  10.  Effect  of  tortuosity  on  the  performance  of  an  anode-supported  SOFC. 


the  porous  layers  increases  the  resistance  to  mass  diffusion  due  to 
increased  diffusion  path  length,  which  results  in  the  reduction  of 
reactant  concentration  at  the  reaction  sites;  as  a  result,  the  contri¬ 
butions  of  activation  and  concentration  overpotentials  to  the  cell 
potential  loss  increases  with  the  increase  of  tortuosity.  Further, 
the  effective  conductivities  of  the  porous  layers  decreases  with  the 
increase  of  tortuosity  resulting  in  the  increase  of  ohmic  contribu¬ 
tion  of  the  porous  layers  to  the  cell  potential  loss,  and  hence  cell 
performance  decreases  with  the  increase  of  tortuosity,  as  shown  in 
Fig.  10. 

The  composition  of  the  electron-conducting  particles  in  the 
reaction  zone  layers  is  an  important  parameter  affecting  the  cell 
performance,  and  its  effect  is  shown  in  Fig.  11.  The  composition  is 
represented  by  the  volume  fraction  of  the  electron-conducting  par¬ 
ticles  (<2>)  in  the  reaction  zone  layers,  which  is  varied  from  0.4  to 
0.6  and  its  effect  on  the  cell  performance  is  illustrated  in  Fig.  11. 
The  other  operating  and  design  conditions  remain  the  same  as  the 
base  case  parameters  given  in  Table  7.  It  is  seen  that  increasing  the 
volume  fraction  of  the  electron-conducting  particles  in  the  reac¬ 
tion  zone  layers  from  0.4  to  0.5  increases  the  performance  of  an 


Fig.  11.  Effect  of  the  volume  fraction  of  the  electron-conducting  particles  in  the 
reaction  zone  layers  on  the  performance  of  an  anode-supported  SOFC. 


anode-supported  SOFC;  however,  further  increase  in  the  volume 
fraction  of  the  electron-conducting  particles  in  the  reaction  zone 
layers  results  in  the  reduction  of  cell  performance.  This  is  due  to  the 
fact  that  the  largest  reactive  surface  area  for  electrochemical  reac¬ 
tions  is  achieved  when  the  dimensions  and  volume  fractions  of  the 
electron-  and  ion-conducting  particles  are  equal  [41  ].  Moreover,  the 
effective  electronic  and  ionic  conductivities  in  the  reaction  zone  lay¬ 
ers  are  functions  of  the  volume  fraction  of  the  electron-conducting 
particles;  increasing  0  beyond  0.5  increases  the  effective  electronic 
conductivity  but  decreases  the  effective  ionic  conductivity  in  the 
reaction  zone  layers  resulting  in  an  increase  in  the  ohmic  overpo¬ 
tential,  and  thereby  decreasing  the  cell  performance. 

6.  Conclusions 

A  general  electrode-electrolyte-assembly  (EEA)  model  has  been 
developed  for  the  performance  characteristics  of  planar  anode- 
supported  solid  oxide  fuel  cells.  One  of  the  novel  features  of  the 
present  model  is  its  treatment  of  the  porous  electrodes  composed  of 
two  distinct  layers,  referred  to  as  the  backing  layer  and  the  reaction 
zone  layer.  The  present  model  is  fuel  flexible  and  includes  methane 
reforming  and  water-gas  shift  reactions  in  the  anode.  The  coupled 
governing  equations  of  species,  charge  and  energy  along  with  the 
constitutive  equations  in  different  layers  of  the  cell  are  solved  using 
finite  volume  method.  The  present  model  can  predict  all  forms  of 
overpotentials  in  the  electrodes  and  the  predicted  anode  concen¬ 
tration  overpotential  is  validated  with  measured  data  available  in 
literature.  An  excellent  agreement  is  obtained  between  the  pre¬ 
dicted  and  measured  concentration  overpotential  with  coefficient 
of  determination  as  high  as  0.989.  It  is  found  that  the  anode  overpo¬ 
tential  is  the  largest  contributor  to  the  cell  potential  loss  at  higher 
current  densities,  followed  by  cathode  and  electrolyte  overpoten¬ 
tials.  Although  the  thicknesses  of  the  cathode  and  the  electrolyte 
in  an  anode-supported  SOFC  are  negligible  relative  to  the  thick¬ 
ness  of  the  anode  but  the  contribution  of  cathode  and  electrolyte 
overpotentials  are  not  negligible.  In  addition,  even  at  low  fuel  uti¬ 
lizations,  the  contribution  of  anode  concentration  overpotential  in 
an  anode-supported  SOFC  becomes  significant  when  there  are  no 
chemical  reactions  in  the  anode.  Flowever,  the  anode  concentra¬ 
tion  overpotential  is  still  three  orders  of  magnitude  smaller  than 
the  anode  ohmic  overpotential.  Reducing  the  operating  tempera¬ 
ture  below  1073  I<  results  in  a  significant  drop  in  the  performance 
of  an  anode-supported  SOFC;  hence  ionic  conductivity  of  the  ion¬ 
conducting  particles  in  the  reaction  zone  layers  and  electrolyte 
need  to  be  enhanced  to  operate  anode-supported  SOFCs  below  1073 
K.  Further,  the  anode  reaction  zone  thickness  of  50  p,m  is  found  to 
be  an  optimal  thickness  for  the  conditions  investigated  in  this  study. 

References 


[1  ]  M.A.  Khaleel,  Z.  Lin,  P.  Singh,  W.  Surdoval,  D.  Collin,  Journal  of  Power  Sources 
130(2004)136-148. 

[2]  S.C.  Singhal,  Solid  State  Ionics  152-153  (2002)  405-410. 

[3]  P.  Aguiar,  C.S.  Adjiman,  N.P.  Brandon,  Journal  of  Power  Sources  138  (2004) 
120-136. 

[4]  H.  Yakabe,  M.  Hishinuma,  M.  Uratani,  Y.  Matsuzaki,  I.  Yasuda,  Journal  of  Power 
Sources  86  (2000)  423-431. 

[5]  J.-W.  Kim,  A.V.  Virkar,  K.-Z.  Fung,  K.  Mehta,  S.C.  Singhal,  Journal  of  the  Electro¬ 
chemical  Society  146  (1999)  69-78. 

[6]  A.V.  Virkar, J.  Chen,  C.W. Tanner, J.-W.  Kim,  Solid  State  Ionics  131  (2000)  189-198. 

[7]  J.H.  Nam,  D.H.  Jeon,  Electrochimica  Acta  51  (2006)  3446-3460. 

[8]  J.  Yuan,  B.  Sunden,  Transactions  of  the  ASME,  Journal  of  Fuel  Cell  Science  and 
Technology  3  (2006)  89-98. 

[9]  J.  Yuan,  B.  Sunden,  Transactions  of  the  ASME,  Journal  of  HeatTransfer  127  (2005) 
1380-1390. 

[10]  T.  Ackmann,  L.G.J.  de  Haart,  W.  Lehnert,  D.  Stolten,  Journal  of  the  Electrochemical 
Society  150  (2003)  A783-A789. 

[11  ]  W.  Lehnert,  J.  Meusinger,  F.  Thom,  Journal  of  Power  Sources  87  (2000)  57-63. 

[12]  J.  Deseure,  Y.  Bultel,  L.  Dessemond,  E.  Siebert,  Electrochimica  Acta  50  (2005) 
2037-2046. 


928 


M.M.  Hussain  et  al.  /  Journal  of  Power  Sources  189  (2009)  916-928 


[13]  X.J.  Chen,  S.H.  Chan,  K.A.  Khor,  Electrochimica  Acta  49  (2004)  1851-1861. 

[14]  S.H.  Chan,  Z.T.  Xia,  Journal  of  the  Electrochemical  Society  148  (2001) 
A388-A394. 

[15]  S.  Sunde,  Journal  of  Electroceramics  5  (2000)  153-182. 

[16]  P.  Costamagna,  P.  Costa,  V.  Antonucci,  Electrochimica  Acta  43  (1998)  375- 
394. 

[17]  R.E.  Williford,  L.A.  Chick,  Surface  Science  547  (2003)  421-437. 

[18]  J.  Fleig,  The  Annual  Review  of  Materials  Research  33  (2003)  361-382. 

[  19]  M.M.  Hussain,  X.  Li,  I.  Dincer,  International  Journal  of  Energy  Research  29  (2005) 
1083-1101. 

[20]  B.A.  Haberman,  J.B.  Young,  International  Journal  of  Heat  and  Mass  Transfer  47 
(2004) 3617-3629. 

[21  ]  B.  Todd,  J.B.  Young,  Journal  of  Power  Sources  110  (2002)  186-200. 

[22]  T.  Kenjo,  S.  Osawa,  K.  Fujikawa,  Journal  of  the  Electrochemical  Society  138  (1991 ) 
349-355. 

[23]  T.  Kenjo,  M.  Nishiya,  Solid  State  Ionics  57  (1992)  295-302. 

[24]  J.  Abel,  A.A.  Kornyshev,  W.  Lehnert,  Journal  of  the  Electrochemical  Society  144 
(1997)4253-4259. 

[25]  D.H.  Jeon,  J.H.  Nam,  C.-J.  Kim,  Journal  of  Power  Sources  139  (2005)  21-29. 

[26]  J.  Newman,  ICE.  Thomas-Alyea,  Electrochemical  Systems,  3rd  ed.,  Wiley,  New 
Jersey,  2004. 

[27]  R.  Bove,  S.  Ubertini,  Journal  of  Power  Sources  159  (2006)  543-559. 

[28]  S.V.  Patankar,  Numerical  Heat  Transfer  and  Fluid  Flow,  McGraw-Hill,  New  York, 
1980. 


[29]  H.K.  Versteeg,  W.  Malalasekera,  An  Introduction  to  Computational  Fluid 
Dynamics:  The  Finite  Volume  Method,  John  Wiley  &  Sons  Inc.,  New  York,  1995. 

[30]  N.R.  Amundson,  E.  Morana,  R.  Sanders,  East-West  Journal  of  Numerical  Math¬ 
ematics  6  (1998)  9-25. 

[31]  J.M.  Stockie,  K.  Promislow,  B.R.  Wetton,  International  Journal  for  Numerical 
Methods  in  Fluids  41  (2003)  577-599. 

[32]  H.  Zhu,  R.J.  Kee,  Journal  of  Power  Sources  117  (2003)  61-74. 

[33]  S.H.  Chan,  C.F.  Low,  O.L.  Ding,  Journal  of  Power  Sources  103  (2002)  188-200. 

[34]  J.R.  Ferguson,  J.M.  Fiard,  R.  Herbin,  Journal  of  Power  Sources  58  (1996)  109-122. 

[35]  M.K.  Alkhateeb,  S.J.  Parulekar,  J.R.  Selman,  S.  Al-Hallaj,  Electrochemical  Society 
Proceedings,  vol.  2005-2007,  2005,  pp.  738-748. 

[36]  N.Q,  Minh,  T.  Takahashi,  Science  and  Technology  of  Ceramic  Fuel  Cells,  Elsevier 
Science  B.V.,  Amsterdam,  The  Netherlands,  1995. 

[37]  M.  Brown,  S.  Primdahl,  M.  Mogensen,  Journal  of  the  Electrochemical  Society 
147(2000)475-485. 

[38]  A.  Bieberle,  L.P.  Meier,  L.J.  Gauckler,  Journal  of  the  Electrochemical  Society  148 
(2001 )  A646-A656. 

[39]  G.O.  Lauvstad,  R.  Tunold,  S.  Sunde,  Journal  of  the  Electrochemical  Society  149 
(2002)  E506-E514. 

[40]  R.E.  Williford,  L.A.  Chick,  G.D.  Maupin,  S.P.  Simner,  Journal  of  the  Electrochem¬ 
ical  Society  150  (2003)  A1067-A1072. 

[41]  P.  Costamagna,  K.  Honegger,  Journal  of  the  Electrochemical  Society  145  (1998) 
3995-4007. 

[42]  I.  Drescher,  W.  Lehnert,  J.  Meusinger,  Electrochimica  Acta  43  (1998)3059-3068. 


